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We study fundamental and vortical solitons in disk-morphed Bose-Einstein condensates (BECs) 
subject to strong confinement along the axial direction. Starting from the three-dimensional (3D) 
Gross-Pitaevskii equation (GPE), we proceed to an effective 2D nonpolynomial Schrodinger equation 
(NPSE) derived by means of the integration over the axial coordinate. Results produced by the 
latter equation are in very good agreement with those obtained from the full 3D GPE, including 
cases when the formal 2D equation with the cubic nonlinearity is unreliable. The 2D NPSE is used 
, to predict the density profiles and dynamical stability of repulsive and attractive BECs with zero 

. and finite topological charge in various planar trapping configurations, including the axisymmetric 

' harmonic confinement and ID periodic potential. In particular, we find a stable dynamical regime 

that was not reported before, viz., periodic splitting and recombination of trapped vortices with 
topological charge 2 or 3 in the self-attractive BEG. 
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INTRODUCTION 



^ A natural setting for the creation of localized states, both fundamental ones and those carrying intrinsic vorticity, 

in Bose-Einstein condensates (BECs) is provided by "pancake" configurations, which are strongly confined in one 
(axial) direction, z, being weakly trapped in the transverse plane, (x,y). In the experiment, "pancakes" were created 
by a superposition of a tight fiat optical trap, formed by a pair of strongly repulsive (blue-detuned) light sheets, and 
^ \ a loose in-plane radial magnetic trap Q • Alternatively, the axial trap may be represented by one cell of a very strong 
. optical lattice (OL) 0- 

Ch ' The pancake configuration may exist with either sign of the intrinsic nonlinearity, self-repulsive or self-attractive 
O , (although, as far as we know, in the latter case experiments were not reported, as yet, in the pancake geometry). In 
I particular, the attractive sign can be induced by means of the Feshbach-resonance technique, as in well-known works 

\ which demonstrated the creation of effectively one-dimensional (ID) solitons in the condensates of ""Li 0] and ^'^Rb 
I ■ [J] atoms, and in other settings (for instance, in experiments with the condensate formed by weakly interacting ^^K 
^ ' atoms (H). Depending on the sign of the nonlinearity, one may expect the establishment of different matter- wave 
■ patterns. In the case of the repulsion, a relatively weak (in comparison with the tight axial trap) in-plane parabolic 
' (alias harmonic) confining potential gives rise to an axisymmetric ground state, that may be described by means of 
. the Thomas- Fermi (TF) approximation . Imparting the angular moment to the confined axisymmetric state creates 
a vortex, which may be unstable against the bending of its axis in the full 3D geometry Q, but is stable in the planar 
form. On the other hand, multiple vortices, with topological charge S > 2, were demonstrated experimentally Q 
and theoretically @ to be unstable against splitting into unitary ones. Nevertheless, it was predicted that multiple 
vortices may be stabilized by imposing an anharmonic in-plane trapping potential [10| . 

In the case of the self-attractive nonlinearity, a challenging problem is to predict, and create in an experiment, 2D 
solitons and solitary vortices (solitons with the embedded vorticity) that would be stable in the pancake setting, despite 
the possibility of the collapse. If the in-plane confinement is provided by the axisymmetric parabolic potential, soliton- 
\^ ' like states are stable provided that their norm (the total number of atoms) falls below a critical value which stipulates 
5^ , the onset of the collapse The critical value strongly increases for vortex solitons, but they are vulnerable, below 
the collapse threshold, to another instability, which tends to break them into a few fragments shaped as fundamental 
solitons, which may later suffer the intrinsic collapse. Stability limits for the trapped vortices in the self-attractive 
BEC have been studied in detail theoretically, both for the nearly planar and full 3D geometries flp |- . including 
the limit case of the cigar-shaped (prolate) trap, opposite to that of the oblate "pancake" [H, Il9l |. In particular, 
it was found that such solitons with multiple charges, S* > 2, are always unstable; as concerns the unitary vortices 
[S = 1), they feature a region of a semi-stability between fully stable and unstable configurations, where the vortex 
periodically splits and recombines [3, |T3| . 

Another general possibility to stabilize both fundamental solitons and various solitary-vortex states in the 2D 
geometry is offered by the use of periodic (rather than parabolic) in-plane potentials, that maybe typically created 
by means of OLs. For 2D periodic potentials, this possibility was first demonstrated in Refs. and [2lj. Then, it 
was shown that a quasi-lD potential, induced by a ID OL (i.e., a periodic potential that depends on a single in-plane 
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coordinate) can also stabilize fundamental solitons [24I. The general case of an anisotropic potential, composed by 
two perpendicular sublattices with unequal strengths, was investigated recently [2^. 

A very accurate model for the description of the mean-field dynamics of dilute condensates is provided by the 3D 
Gross-Pitaevskii equation (GPE) [6|. An effective 2D equation for the "pancake" geometry should be derived from 
the 3D GPE by means of an adequate reduction procedure. It was developed, using the ansatz which assumes the 
factorization of the wave function in the axial direction and perpendicular (pancake's) plane, in Ref. [23| (see also Ref. 
[2^). The result is a nonpolynomial nonlinear Schrodingcr equation(NPSE) (it takes a different nonpolynomial form 
in the cigar-shaped setting, when the dimension is reduced from 3 to 1 [l^-jlgl). The objective of the present work is 
to use the 2D NPSE for the systematic analysis of effectively 2D localized states in both cases of the self-repulsive and 
self-attractive nonlinearity. To this end, we briefly recapitulate the derivation of the effective 2D equation in Section 
II, where the respective TF approximation is considered too. In Section III we focus on the model with the repulsive 
nonlinearity. Using the 2D NPSE, we construct solution families for both the ground state and vortices trapped in 
the axisymmetric in-plane parabolic potential. The comparison with direct numerical solutions of the underlying 3D 
GPE demonstrates that the effective equation, unlike the formal 2D version of the GPE, with the cubic nonlinearity, 
provides for virtually exact shapes of both the ground and vortex states. In Section IV, we consider the model with the 
self-attraction. First, we construct families of fundamental and vortical states trapped in the axisymmetric parabolic 
potential. In that case too, we conclude that the respective version of the 2D NPSE provides for a virtually exact 
prediction of the shape of the states. The 2D NPSE reproduces earlier known results for the stability of trapped 
states with vorticity S ~ 1, and generates a new stable regime for S = 2 and 3, in the form of periodic splitting of the 
double or triple vortex into a pair or triplet of unitary vortices followed by their recombination back into the single 
double or triple vortex. Then, we consider the quasi-lD periodic potential, for which families of fundamental solitons 
are found as numerical solutions of the 2D NPSE, and also by means of the variational approximation (VA) applied 
directly to the underlying 3D GPE. 



II. THE EFFECTIVE 2D EQUATION 



A. The general case 



We consider the condensate formed by a large number N of bosonic atoms of mass m tightly trapped in axial 
direction z by the harmonic potential, {m/2)ujlz'^, and loosely confined in transverse plane (a;, y) by a generic potential, 
W{x,y). The scaled 3D GPE, which governs the mean-field dynamics of the dilute potential at zero temperature, is 



■9 , 



(1) 



where V{x, y, z) = W{x, y) -|-z^/2, g = 2Nas/az is the interaction strength, with the inter-atomic scattering length 
and Qz = ^JhJimjjJz) the scale of the axial trapping, and 3D wave function ^'(r, t) is subject to normalization 

J\^{r,t)\Uh = l. (2) 

Following Ref. [li], the solution to Eq. (P) can be simplified by setting 

*(r,t) = $(x,y,t) S(z;77), (3) 
where is an arbitrary transverse wave function, and ^(z,?]) is its normalized axial counterpart, 

S(z;^) = (7rry2)-i/%xp[-zV(2ry2)], / E^{z)dz^l. (4) 

J —00 

With regard to underlying condition the factorization ansatz based on Eqs. ^ and ^ yield the normalization 
condition for the 2D wave function, J J \^(x,y,t)\'^ dxdy = 1. In these expressions, axial width r] may depend on 
transverse coordinates {x, y) and time t. 

The substitution of factorized expression ([3]) in Eq. ([1]) and subsequent averaging in axial coordinate z lead to the 
following 2D NPSE [H], 



■9 ^ 



2 77 4 \?7^ 



(5) 



3 



2 



1 1 1 1 1 1 




exact 






strong- attractive 






strong-repulsive 




j 
/ 


1 


1 



100 



-50 



50 



100 



FIG. 1: (Color online). Axial width rj of the pancake-shaped condensate as a function of the interaction-strength x density 
product, 7|4>p. The solid curve: the full dependence following from Eq. (|8|; dashed and dotted curves: asymptotic expressions 
(|12l) . in the cases of repulsion, 7 > 0, and attraction, 7 < (except for the narrow vicinity of 7|"I'P = 0, the asymptotic curves 
are virtually identical to the full one, in the regions of 7|4'P > and 71$^ < 0, respectively). 



where the scaled interaction strength is 

7 = 2V2TTNas/a:, = V^ng, 
and the axial width is determined by an algebraic equation. 

Exact solutions to Eq. ([7]) are given by the Cardano formula, 



where the upper and lower signs correspond, respectively, to 7 > and 7 < 0, and 

1 /3 

A = (3/2)1/^ (97'!$!^ + \/3V256 + 27-/^\<P\»^ . 
The dependence of 77 on product 7|<I>P, as produced by Eqs. ^ and (O, is displayed in Fig. [TJ 
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B. Low- and high-density limits 



In the low-density limit, |7||$P <C 1, Eq. ([7]) yields a simple solution. 



1+7|$|V4. 



(10) 



The substitution of this approximation in Eq. ([5|) reduces it to the nonlinear Schrodinger equation (NPSE) with the 
cubic-quintic (CQ) nonlinearity, in which the quintic term always corresponds to the effective self-attraction: 
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-ivi + W(x,j/)+7|$p-^7'l$r 



(11) 



The derivation of the effective ID equation from the 3D GPE in the low-density limit also leads to the NPSE with 
the CQ nonlinearity of the same type [131 . Note that the self- focusing quintic term in the 2D NPSE leads to the 
supercritical collapse (the one with zero threshold) . If the 2D equation includes both quintic and cubic self-attraction 
terms, then a periodic potential, accounted for by W(x. y) in Eq. ()lip . can stabilize various types of solitons and 
localized vortices against the (supercritical) collapse (28|. However, the 2D equation with the quintic nonlinearity 
only, see Eq. ([T^ below, cannot give rise to stable solitons, irrespective of the presence of the periodic potential [28j . 



4 



In the high-density Umit. |7||3'|^ ^ 1, asymptotic expressions for the axial width which follow from Eq. ([7]) are 
different in the cases of repulsion (7 > 0) and attraction (7 < 0), 

^S-(7l</>n^/^e^ = -(7l<fr)". (12) 

The substitution of these approximations in Eq. ([5]) leads to two different 2D NPSEs: in the case of the repulsion, it 
is 



'at''- 



-ivi + M/(x,y) + j7'/3|$|4/3 



(13) 



which coincides with the equation of the hydrodynamic type for the degenerate fermion gas in the weak-coupling 
regime [2^, and in the case of the attraction, the effective equation is the NPSE with the quintic term only: 



1 fd^ 



(14) 



C. The Thomas-Fermi approximation 

The substitution of ^{x,y,t) ~ (j){x,y) e~''^*, with chemical potential p., casts Eq. ([5]) in the stationary form, 
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= M(a;,y)(^, 



(15) 



where /i(x, y) = p. ~ W{x, y) is the local chemical potential. In the Thomas-Fermi (TF) approximation, that may be 
relevant in the case of the repulsive interactions, i.e., 7 > 0, one neglects the spatial derivatives in Eq. (jlSp . and, 
taking into regard Eq. ([8]), obtains the following expression. 



I0(a;,y)|' 
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/y^||^Kx,2/)]' + i^+/.(x,2/") 



(16) 



which is valid as long as it yields |0(a:,y)|2 > 0. In the region where expression (|16p is negative, it is replaced 
by \(t>{x,y)\'^ = 0. In the case of |/i (x, y) — 1/2| <C 1/2, the TF approximation corresponds to the low-density 
limit PU)) . with \(j){x,y)\^ « 7"-^ [^(a;, y) — 1/2] . The high-density limit corresponds to large which yields 
|</.(x,y)p«7-l[(4/5)/i(x,y)]^/^ 

III. THE GROUND STATE AND VORTICES IN THE 2D PARABOLIC TRAP, WITH THE 

REPULSIVE NONLINEARITY 



We start the analysis of trapped states, predicted by the 2D NPSE, with the case of the repulsion, assuming that 
the potential inside the "pancake" is also a parabolic trap (which is much weaker than the underlying axial trapping 
potential, /2), i.e.. 



W(x,2/) = (AV2)(x2+y2) 



(17) 



in Eqs. ([5]) and pS]) . with small A. For comparison of localized states predicted by stationary NPSE (jlSp . and by the 
full 3D GPE, we introduce the radial and axial probability densities, defined in terms of the 3D wave function, which 
is subject to normalization 



P3d(?', t) 



+ CXD 



|*(r,t)|2 dz, p3B{z,t) = 27r / \^{r,t)\\dr, 



(18) 



where r = ^/ x'^ + y^- In the framework of the 2D description, the factorized ansatz based on Eqs. ([3]) and ([4]) yields 



P2D(r,i) = |$(r,t)P, p2u{z,t)^2V^ 



Vir) 



(19) 
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FIG. 2: (Color online). Radial (a) and axial (b) probability densities, p(r) and p{z), for the self-repulsive condensate in the 
ground state trapped in parabolic potential (|17p with A = 0.1. Solid, dashed, dotted, and dotted-dashed curves display the 
results produced, respectively, by the 3D GPE, 2D NPSE, the ordinary 2D equation with the cubic nonlinearity, and TF 
approximation to the 3D GPE. Values of interaction strength g = 2Nas/az are indicated in the panels. 



In Fig. [51 we display typical examples of the radial and axial densities for the repulsive condensate in the ground 
state, fixing A = 0.1. The respective 3D and 2D stationary wave functions were obtained by means of the imaginary- 
time integration, using the finite-difference Crank-Nicolson predictor-corrector code, as elaborated in Ref. (30(|. The 
figure demonstrates that the 2D NPSE provides almost exact results, in the comparison to the full 3D equation, except 
for some deviation of the axial density at very large values of the interaction strength, g = 200 and 1000. It is worthy 
to note that the NPSE provides for a much better accuracy than the formal 2D GPE with the cubic nonlinearity. 

In the same setting, the 3D wave functions for vortex states can be sought for as ^(r, = ip{r,z,t) e*"^^, where r 
and 9 are the polar coordinates in the (x, y) plane, S is integer vorticity, and function ip obeys equation 



ld_ 

r dr 



+A^r^ + z-' + 47rg|V'| 



+ 



(20) 



In the case of 2D NPSE ([5]), the vortex corresponds to <E'(r, 6, t) 

.3 , 1 



■dt 



ld_ 

r dr 



r,t)e'^^, with Eq. replaced by 



+X\^ + 27^ + 

T] 2 



V 



(21) 



Relevant solutions to both equations () and ([51]) must feature the usual asymptotic form, (■;/;, 0) ~ r'^', at r ^ 0. 

The radial and axial profiles of the vortex states with 5 = 1 are displayed in Fig. O for the same value of A = 0.1 
as in Fig. [21 To generate these results, Eqs. (j^Ol) and (|5T|) were also solved in the imaginary time by means of the 
above-mentioned finite-difference Crank-Nicolson algorithm [s^. We again conclude that, as well as in the case of the 
ground state, the 2D NPSE produces results which, in most cases, virtually coincide with those obtained from the 
3D GPE (except for some deviation of the axial profiles at very large g), while the formal 2D GPE with the cubic 
nonlinear term is essentially less accurate. 
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For higher values of S (up to S = 3), the radial and axial density profiles, as predicted by the 2D NPSE, are 
displayed, along with their counterparts for the ground state, in Fig. 31 at fixed values of the parameters, A = 0.1 
and g = 200. As expected, the growth of S leads to the increase of the size of the vortex' core. On the other hand, 
the axial density profile, p{z), is not strongly affected by the change of S. 

We have checked that, as well as in the case of full 3D GPE with the repulsive nonlinearity and its 2D counterpart, 
multiple trapped vortices are unstable against splitting into solitary vortices within the framework of the 2D NPSE. 
This finding is not surprising, as the latter equation is derived from the 3D GPE. hence it should inherit basic 
properties of solutions of the 3D equation. 



IV. SOLITONS AND SOLITARY VORTICES IN THE CASE OF THE ATTRACTIVE NONLINEARITY 



A. The axisymmetric parabolic trap 



In the case of the attractive nonlinearity, 5 < 0, localized states and vortices supported by the in-plane potential 
- in particular, axisymmetric parabolic potential p?)) - may be considered as solitons In Fig. [5] we display 

the radial and axial profiles of the solitons found from the NPSE with potential ([T7|) . for A = 0.1 and g = —1. It can 
be checked that, for sufhcicntly small values of \g\, at which the collapse does not take place, the full 3D GPE yields 
virtually the same profiles. 

Under normalization conditions ([2]) and parabolic confinement (|17p . the critical value, gc, of the interaction strength 
for the onset of collapse can be found as that at which the stationary soliton solution disappears. These values, as 
found from the set of the three models (3D GPE, 2D NPSE, and the formal 2D GPE, with A = 0.1), are collected in 
Table 1, for different values of vorticity S. Naturally, gc quickly increases with S [15|, and we once again conclude 
that the NPSE provides for high accuracy, especially in comparison with the large error in the predictions produced 
by the formal 2D equation with the attractive cubic nonlinearity. 
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( 5c)3D-GPE 


(5c)2D-NPSE 


(5c)2D-GPE 





-1.9 


-1.9 


-2.3 


1 


-7.3 


-7.2 


-8.9 


2 


-12.4 


-12.3 


-15.8 


3 


-16.9 


-16.8 


-21.6 



Table 1. Critical strength gc for the onset of the collapse in the attractive condensate with vorticity S, confined by 
in-plane potential (fTT]) with A = 0.1. The table compares the predictions produced by the full 3D GPE, the effective 
2D equation (NPSE), and the formal 2D equation (GPE) with the cubic nonlinearity. 




FIG. 3: 



(Color online). The same as in Fig. [21 but for vortex states with S = 



1. 
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FIG. 4: (Color online). The upper and lower panels display, respectively, the radial and axial probability densities, p(r) and 
p{z), for the ground state and higher-order vortices, as found from the 2D NPSE with A = 0.1 and g = 2Nas/a^ — 200. The 
solid, dotted, dashed, and dotted-dashed curves pertain, respectively, to S — 0, S = 1, S — 2, and S = 3 (all vortices with 
S > 2 are unstable against splitting into unitary vortices). 
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FIG. 5: (Color online). The same as in Fig. |4] (again, with A — 0.1), but for the self-attractive condensate, with nonlinearity 
strength g = 2Naa/a^ = — 1. 



An important issue is the stability of solitary vortices with vorticity S in the axisymmetric harmonic trap. We 
tackled this problem by means of direct simulations of the time-dependent 2D NPSE, Eq. ([5]), using a finite-difference 
real-time Crank-Nicolson algorithm [sot . The initial conditions were taken as 



$(x, y, t = 0) = C (x + iy)^ cxp 



r 

(52 



(22) 



where C is a normalization constant, and 5^1 accounts for possible breakup of the axial symmetry of the vortex, i.e., 
an initial perturbation applied to the vortex. Note that, with (5 = 1, Eq. ([22]) gives the exact quantum- mechanical 
wave function of the stationary vortex configuration (for 5 = 0). 

We report results of the simulations for 5 = 1.1 and trap's strength A = 0.1. Wc have concluded that the trapped 
vortex with = 1 is stable for g > —1.3: in this region, the initial perturbation leads only to intrinsic periodic 
oscillations of the vortex, which maintains its integrity. In agreement with previously published results fTsj-jlTj. the 
value of the nonlinearity strength at the stability border, ^stab = —1.3, is much lower than at the collapse threshold, 

gc = -7.2. 

Outside the stability domain, i.e., at —g > 1.3, wc observed splitting of the trapped vortex into two solitary 
fragments, which would rotate and recombinc back into the vortex. In the interval 1.3 < —g < 1.4 the splitting- 
recombination cycle repeats itself quasi-pcriodically, as shown in Fig. [51 where we have plotted the planar density 
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FIG. 6: Density profiles \${x,y,t)\'^ for tfie solitary vortex with S = 1 and g — —1.4, obtained by solving the time-dependent 
2D NPSE, Eq. with A = 0.1. The four snapshots correspond to t = (a), t = 56 (b), t ^ 88 (c), t = 120 (d). Time is given 
in units of l/uj^- The figure displays the quasiperiodic process of the alternating splittings and recombinations of the vortex 
ending with the collapse, see the text. 



profile, at four instants of time correspon ding to different stages of the cycle. Note that a similar 

dynamical regime was previously demonstrated in Refs. [l3l. [l5|. in the framework of the ordinary 2D equation with 
the cubic self-attractive term. If the self-attraction is stronger, viz., 1.4 < —g < 1.5, the two fragments suffer an 
intrinsic collapse after two cycles (in the case shown in Fig. [6]thc collapse takes place at t ~ 125). For a still stronger 
self-attraction, —g > 1.5, the two rotating fragments do not recombine at all, quickly undergoing the intrinsic collapse. 

We have also investigated the stability of the solitary vortex with 5 = 2. In agreement with Ref. [lB|, it was 
found that the double vortex can never be stable against splitting into a pair of unitary ones. Nevertheless, a stable 
dynamical regime, which has not been reported in previous studies, was found in this case: the simulations of Eq. 
(O, with initial conditions (|^^ corresponding to S* = 2, demonstrate that, at small values of the interaction strength, 
—g < 0.5, the double vortex periodically features incomplete separation into two vortices, which manifests itself in the 
splitting of the zero-amplitude point (the center of the vortex) into a pair of local minima, separated by a relatively 
small distance [see Figs. [7l^b,d)], which is followed by recombination of the pair back into a single point, see a set of 
snapshots in Fig. [T] 

At larger values of the interaction strength, —g > 0.5, the switchings between the states with one and two local 
minima is a transient regime, which is followed by permanent splitting of the double vortex into a pair of unitary 
vortices rotating around point a; = y = 0, as shown in Fig. [Sj Finally, at —g > 1.8 the evolution of the double vortex 
ends up with the collapse. 

For the sake of the completeness of the analysis, we have also analyzed the dynamics of the triple vortex, with 
5 = 3. It features a trend to the splitting into a state with three local minima of the amplitude. At small enough 
values of \g\, the splitting is periodically followed by the recombination into the original triple vortex, while at larger 
values of \g\ it permanently splits into a set of three rotating unitary vortices. 
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(a) 



(b) 





(c) 



(d) 





FIG. 7: The same as in Fig. [G] but for the solitary vortex with S — 2 and g — —0.3. The snapshots pertain to t = (a), 
t = 16 (b), t = 32 (c), t = 48 (d). The figure displays one-and-a-half cycles of the periodic process (of period r ~ 32) of the 
alternating splittings and recombinations of the vortex, see the text. 



B. The periodic in-plane potential 



As discussed in Introduction, a problem of considerable interest is the stabilization of fundamental (non-vortical) 
2D solitons in the attractive condensate loaded into a quasi-lD periodic potential, corresponding to a ID OL (optical 
lattice); the same problem has another physical realization in terms of nonlinear optics [l^]. In this section, we 
address the problem within the framework of NPSE ^ with potential W{x, y) = — Vq cos {2kx). Accordingly, the full 
potential in 3D GPE (P) is V{x, y, z) = -Vq cos {2kx) + /2. 

First, the corresponding stationary solutions to the full three-dimensional GPE, Eq. ([1]), are looked for as 'I'(r, t) ~ 
ipir) exp(— i/it), with chemical potential fi and real function ip which obeys equation 



-iv^ - V"oCOs(2A:x) + ^z^ + 2TTg\i'\ 



ip = flip. 



This equation can be obtained by minimizing the energy functional, 

E= [ (fr iP* 



-iv^ - Vq cos {2kx) + ^z'^tp + TTg\ip\ 



subject to constraint ([J), which results in the following expression for the chemical potential, 



(23) 



(24) 



(25) 



To predict solitons by means of the VA (variational approximation), we use the normalized 3D Gaussian ansatz 
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FIG. 8: The same as in Fig. □ but for g = -0.8. The snapshots pertain to t = (a), f = 44 (b), t = 136 (c), t = 168 (d). The 
figure displays the permanent sphtting of the double vortex into a pair of rotating unitary vortices, see the text. 
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with widths ci, (T2 and a^. Inserting the ansatz into Eq. (|24[) . one obtains 
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(27) 



and Eq. (j25p yields, in the same approximation, fi ~ E + g/ (v27r(7i(T20'3) ■ Aiming to predict the ground state 
in the framework of the VA, we look for values of tri, (T2, and CT3 that minimize the energy, imposing conditions 
dE/dai^2.3 ~ 0. In this way, we derive the following variational equations. 



4VbfcVi exp(-fcV^) 
J_ ^ 9 1_ . 

0'2 \/27rO'l(T3 
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CTl 



2tt o-2<J3 
9 era 



2tt (Jicr2 



(28) 



A solution corresponds to a true minimum of the energy provided that the respective Hessian, d^E/daidaj, is positive 
definite. 

Analysis of the variational equations demonstrates that they may yield the ground-state solution if the OL is strong 
enough, viz., Vq > e^k^/lG. Under this condition, the solutions exists for < 5 < 0, where critical strength 9c 
depends on Vq and k. Notice that, restricting the VA to the 2D setting, i.e., fixing 0-3 ~ 1, cf. Eq. one will instead 
find stable solutions for 9c < 9 < 9s, with 9c — — V27r, and gs — — ^27r [1 — l6Vo/{e'^k^)] for Vq < e^fc^/16, while 
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FIG. 9: (Color online). In the upper panel, the dashed curves depict widths ai, a"2, era of the soliton, in the model with the 
quasi-lD in-plane potential, versus self-attraction strength g, as predicted by the VA (variational approximation) for Vo = k = 1. 
The solid lines with dots show the widths as found from a numerical solutions of the 2D NPSE. The lower panel displays critical 
interaction strength for the onset of the collapse versus amplitude Vo of the periodic potential, as predicted by the VA, for 
two different values of the wavenumber of the periodic potential, k = 1 and k = 2. 



= for Vo > [24|. The expression for predicts the minimum absolute value of the self-attraction strength, 

above which the state cannot remain localized in the periodic potential (in fact, the vanishing of gs for Vq > e^fc^/1 is 
an artifact of the VA 
never vanishes 2^, 



22, 



22 , 31 1 ; in the numerical solution reported for the 2D GPE with the self-attractive cubic term. 



31|). 



Numerically found solutions of Eqs. are plotted versus g in the upper panel of Fig. [SI for Vq = 1 and k = 1 
(dashed lines). We include also the corresponding widths obtained by solving the 2D NPSE. The lower panel of Fig. 
HI shows critical interaction strength gc versus amplitude Vq of the periodic potential. At \g\ > \gc\, solutions for 
solitons do not exist anymore. 

Comparing the predictions of the VA with results produced by the 2D NPSE, we conclude that, in the case of 
weak nonlinearity, the VA based on the Gaussian ansatz is wrong, which is not surprising: in that case, solitons 
are not localized like Gaussians, but rather extend over several cells of the periodic potential [13, [3l|. In particular, 
similar to what was found in the 2D model with the cubic nonlinearity, and contrary to the prediction of the VA, the 
above-mentioned critical value gs, below which (at \g\ < \gc\) the dclocalization of the 2D soliton takes place, never 
drops to zero, as shown in Table 2. On the other hand, critical strength gc for the onset of the collapse, as predicted 
by the VA, is indeed close to what is found from the numerical solution of the 2D NPSE, which is explained by the 
fact that the collapse occurs in well localized configurations. Direct simulations of the time-dependent 2D NPSE 
corroborate a natural expectations that the fundamental soliton is always stable in its existence range, gc < g < gs- 



Vo 


9c 


9s 


1 


-1.12 


-0.91 


2 


-0.92 


-0.43 


3 


-0.81 


-0.25 



Table 2. Critical strengths gc and g^. between which the fundamental soliton exists (and is stable) in the model with 
the quasi-lD periodic potential with amplitude Vq and wavenumber fc = 1, as found from the numerical solution of 
the 2D NPSE. 

Finally, typical profiles of the soliton in the present version of the model, as found from a numerical solution of Eq. 
are shown in Fig. [TOl by means of the density integrated in the free (y) direction, p{x) = ip"^ {x, y) dy. The 
figure corroborates that the profile indeed features a trend to dclocalization with the decrease of \g\. 
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FIG. 10: (Color online). Density profiles of the fundamental solitons stabilized by the quasi-lD periodic in-plane potential are 
shown for three values of self-attraction strength g, while the amplitude and wavenumber of the potential are fixed, Vb = —2 
and k = 1. The density, p{x), is integrated in the free direction, y. 



In this work, we have reported results of the analysis of localized states in nearly-2D ("pancake" -shaped) BEC, with 
both signs of the intrinsic nonlinearity, repulsive and attractive. The localized states represent both the ground state, 
with zero vorticity, and excited states in the form of localized vortices. The main objective of the work was to apply 
the nonpolynomial nonlinear Schrodinger equation (NPSE) derived from the full 3D GPE (Gross-Pitaevskii equation), 
to the analysis of the pancake-morphed localized states. First of all, it has been demonstrated that, in the case of the 
self-repulsive nonlinearity and relatively weak parabolic (harmonic) in-plane axisymmetric trapping potential, the 2D 
NPSE predicts both the radial and axial profiles of the ground states and vortices, with the topological charge up to 
5 = 3, in a virtually exact form, if compared with numerical solutions of the underlying 3D GPE. On the contrary to 
that, the formal application of the ordinary 2D equation with the cubic nonlinearity, or Thomas- Fermi approximation 
derived form the 3D GPE, yields a large error. Similar conclusions have been made for the fundamental and vortical 
solitons in the model with the same axisymmetric trapping potential and self-attractive nonlinearity. As concerns 
the stability, a new result was reported for the trapped vortices with 5 = 2 and 3 in the self-attractive potential: at 
sufficiently small values of the nonlinearity strength, they feature a stable dynamical regime in the form of periodic 
splittings (respectively, into a pair or triplet of unitary vortices) and recombinations. 

Finally, fundamental solitons were studied in the model combining the ID in-plane periodic trapping potential, which 
does not depend on one of the planar coordinates, and the self-attractive nonlinearity. In this case, we compared 
predictions of the VA (variational approximation), derived for the 3D GPE, and numerical solutions of the 2D GPE. 
As a result, it has been concluded (as might be expected) that the VA yields reasonable results, provided that the 
nonlinearity is not too weak, in which case the underlying assumption of the localization of the fundamental soliton 
around one cell of the periodic potential is not valid. 
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